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We consider modifications of the standard nondinear QCD evolution in an attempt to 
account for some of the missing ingredients discussed recently, such as correlations, discrete- 
ness in gluon emission and Pomeron loops. The evolution is numerically performed using the 
Balitsky-Kovchegov equation on individual configurations defined by a given initial value of 
the saturation scale, for reduced rapidities y = (a s N c /ir) Y < 10. We consider the effects of 
averaging over configurations as a way to implement correlations, using three types of Gaus- 
sian averaging around a mean saturation scale. Further, we heuristically mimic discreteness 
in gluon emission by considering a modified evolution in which the tails of the gluon distri- 
butions are cut-off. The approach to scaling and the behavior of the saturation scale with 
rapidity in these modified evolutions are studied and compared with the standard mean-field 
results. For the large but finite values of rapidity explored, no strong quantitative difference 
in scaling for transverse momenta around the saturation scale is observed. At larger trans- 
verse momenta, the influence of the modifications in the evolution seems most noticeable in 
the first steps of the evolution. No influence on the rapidity behavior of the saturation scale 
due to the averaging procedure is found. In the cut-off evolution the rapidity evolution of 
the saturation scale is slowed down and strongly depends on the value of the cut-off. Our 
results stress the need to go beyond simple modifications of evolution by developing proper 
theoretical tools that implement such recently discussed ingredients. 
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1 Introduction 



The B- JIMWLK 1 equations, derived over the last decade as a result of the concerted 
effort of several groups [1—15] , build upon the original ideas on gluon saturation set 
out in [16] to provide the contemporary description of the evolution of QCD scattering 
amplitudes at high energy. Unfortunately, these are complicated equations and their 
complete solution is unknown. A combination of numerical [17-25] and analytical [26- 
32] studies has established the asymptotic properties of the B- JIMWLK equations. 
Most of these studies have focused on the mean field limit of the equations, where the B- 
JIMWLK set reduces to a single closed equation [12,33] - the Balitsky-Kovchegov (BK) 
equation. Nevertheless, the results obtained for the BK equation are known to deviate 
at the most ~ 10% from those obtained numerically for the full B- JIMWLK [25]. 

In obtaining this set of equations, the gluonic density in the target was assumed 
to be large, and the projectile was taken as a dilute object. The strict validity of the 
resulting evolution scheme is, therefore, restricted to this physical domain. Further, 
it has become clear that the physical region of interest - i.e. that within reach of 
accelerators that are currently operating and those that will be operational in the near 
future - lies within the pre-asymptotic region of the evolution and that in this region the 
evolution is dominated by the, as yet poorly understood, initial conditions [24,34,35]. 

In the last two years, there has been a spurt of activity in this domain triggered 
by the crucial observation that important effects for the evolution are absent from the 
currently available evolution scheme [36-39]. These new contributions stem from the 
understanding of the importance of gluon fluctuations and have been variably referred 
to as 'Pomeron loops' [39-46], 'fluctuations' [38] or 'wave function saturation effects'. 
A particularly important development has been the finding of a duality transformation 
connecting the high and low density regimes [45,46]. This duality places strict con- 
straints on the form of the evolution kernel in the intermediate region. Nonetheless, 
the evolution in this intermediate domain is not completely understood. 

Until now all analysis of the effect of these new ingredients have been addressed, 
both analytically and numerically, through the analogy of high-energy QCD evolution 
to a reaction-diffusion process [38]. These discussions are based on the study of the 
stochastic Fisher-Kolmogorov-Petrovsky-Piscounov (sFKPP) equation [39,42,47-50]. 
1 Balitsky-Jalilian-Marian Iancu McLcrran Wcigcrt-Lconidov-Kovner. 
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The main conclusions extracted from these studies are: i) at small rapidities, mean 
field evolution through the BK equation remains approximately valid, but the rapidity 
evolution of the saturation scale is slowed down and the domain of geometric scaling is 
restricted, in comparison with the results of pure BK evolution; ii) at large rapidities, 
geometric scaling is destroyed, being replaced by a new type of scaling, sometimes 
referred to as diffusive [51]. 

At this stage, and while efforts towards the derivation of the complete evolution 
take place, it appears possible to attempt to implement some of the effects as effective 
modifications of the BK equation, and to examine whether the results agree quali- 
tatively or not with our current knowledge and expectations from the results of the, 
previously mentioned, full high-energy QCD evolution. In this paper, we assess the 
impact of some of the proposed new contributions. In particular, we address the effect 
of correlations by introducing an averaging procedure for the initial conditions for the 
evolution. Further, we model possible effects due to the discreteness in gluon emission 
by considering the BK evolution of gluon distributions with cut-off tails. 

The paper is organized as follows. In Section 2 we outline the numerical method 
used to implement the BK equation, and the procedures adopted to introduce modi- 
fications in the evolution intended to mimic the effects of correlations and gluon dis- 
creteness. Our results are presented in Section 3 and discussed in Section 4. 

2 Numerical method 

Non-linear QCD evolution is implemented using the BK equation [12,33] in momentum 
space in the local approximation (i.e. neglecting the impact parameter dependence). 
This approximation should be valid for a large nucleus in the region far from its pe- 
riphery. The evolved object is the unintegrated gluon distribution (f>(k, b) = 4>{k) with 
k the transverse momentum, b the impact parameter and where the dependence on 
rapidity Y is implicitly understood. The unintegrated gluon distribution is related 
to the scattering probability of a qq dipole of transverse size r on a hadronic target, 
N(r,b) = N(r), through 




(1) 
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For the azimuthally independent piece of the solution, which gives the dominant con- 
tribution at large rapidities, the BK equation reads 

8(f>(k) f dk'~ 



dy J k' 2 



k a (j){k') - k 2 (f)(k) k 2 (f)(k) 
\k' 2 -k 2 \ + y/Ak' A + k 4 



m, (2) 



with the reduced rapidity y = (a s N c /7r) Y. In this equation a s is fixed. For a s = 0.2, 
the maximum y = 10 considered in this paper corresponds to large physical rapidities 
Y ~ 50. 

The equation is solved using a 4th-order Runge-Kutta algorithm with step Ay = 
0.025. The integral is evaluated in the domain —15 < Ink < 35 using a Gauss- 
Chebyshev quadrature with 400 points. Henceforth units of transverse momenta will 
be GeV, with the units corresponding to other quantities to be read from the respective 
equations. By varying the integration limits, the number of points, the step Ay, and by 
comparison with previous algorithms [17,19,23], the accuracy of the computations can 
be estimated to be better than 1 % in a wide region excluding one order of magnitude 
from the limits of the domain. It is usually much better than that. The size of this 
domain allows us to safely study the evolution up to y — 10. 



2.1 Initial conditions 

We will use Eq. (J2J) to evolve individual configurations at initial rapidity y = 0, 
characterized by some functional form and by a given value of the saturation scale 
Q s0 = Q s (y = 0). The first functional form (GBW) essayed is motivated by the 
Golec-Biernat-Wusthoff model [52], 



1 f k 2 

2 Ei IqL 



W(*) = -- Ei ) , (3) 



with Ei the exponential-integral function. The second functional form (MF) reads 

<WAO = 7£ + r(o,o + m£, e=(r?) , (4) 



with 7£ the Euler constant and T the incomplete gamma function. This form is mo- 
tivated by phenomenological studies of geometric scaling in 7*-nucleon and nucleus 
collisions [34]. It behaves oc £ (ln£) for k ^> (^C)Qso- The parameter 5 in this func- 
tional form acts as one minus the anomalous dimension governing the behavior of <fi 



4 



for large k. We take 2 5—1. 

Values of the saturation scale in the region 1CT 2 < Q 2 s0 < 10 2 have been explored 
for both initial conditions. This region defines our averaging domain and has been 
sampled in 284 points which are roughly equidistant in logarithmic scale. 



2.2 Averaging procedure and cut-off 

After evolving up to a given rapidity y, we compute the average of <p(k) as proposed 
in [53] i.e. 



■ C-2 clQlW(Q* )<p( k) 
Ao- dQl W(Q%) 

is performed. Two different averaging procedures have been considered: the linear 
averaging [53] 

(<&-<<&>) 2 " 



W 1 (Qi )=exp ^ 
and a logarithmic averaging procedure given by 



(lnQ 2 -ln(Q 2 ))^ 



(6) 



(7) 



QsO 

with (Q 2 s0 ) = 1. Values of A, A^ = 0.01, 0.1, 1, 10 and 100 have been explored. Let 

us note that for both for the linear and the logarithmic case, A, A; = 0.01 is almost 

equivalent to a 5-function, while for A, A/ = 100 roughly 50 % of the normalization 

of the function is contained in the region we study, 10~ 2 < Q 2 s0 < 10 2 . Finally, W 2 

with Ai = 100 is essentially a flat function, with a variation ~ 20 % in the considered 

region - thus it is our model case for a wide averaging. Let us also indicate that 

although the practical relevance of such large values of A, A; may be questionable, our 

aim in this paper is to examine whether such a modification of standard BK evolution 

has any effect or not - thus our attitude of taking extreme values of the parameters 

characterizing the width of the weight functions for the averaging. 

2 Initial values S > S c , S c = 0.627..., are known [30] to develop a wave front with 5 = S c . On 
the other hand, initial values S < S c produce wave fronts which preserve this initial value of S, a 
behavior which has been numerically verified in [48] and apparently contradicts previous numerical 
studies [24]. We have checked that the reason for this apparent contradiction lies in the different 
regions of transverse momentum studied in both references. The wave front is characterized by the 
initial value of 6 for very large values of k studied in [48] but not considered in [24] . 
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Finally, motivated by the results obtained in the framework of the sFKPP equa- 
tion [38,39,42,47-50] which indicate that at large rapidities the dispersion of a set 
of initial conditions becomes linearly dependent on the rapidity, we will examine the 
effects of a logarithmic averaging procedure W2(Q 2 ) but with a rapidity dependent 
value of the dispersion A; = y. Again, the choice of this value is motivated not by 
phenomenological considerations but by the requirements of visibility of the proposed 
modification in our limited space of sampled initial configurations, see below. 

On the other hand, and following the suggestion about the discreteness of the 
evolution in [38] , see also [47, 48] , we have introduced what we call a cut-off evolution: 
at every step in the evolution in Eq. (J2J), values of 4>{k) < k in the r.h.s. of this equation 
have been set to 0: 

d(f>(k) r dk' 2 lk' 2 <j){k')6[<j){k') - «] - k 2 <j){k)6[(j){k) - k] k 2 <p{k)e\(t){k) - «]" 
dy ~ J ~W \k' 2 - W\ + V4fc /4 + k A 

- <t>(k) 2 0[<t>(k) - «], (8) 

with 9[x] the step function. Values of the cut-off k = 0.002, 0.01 and 0.05 have been 
used, for evolution starting from a given initial condition with Q s0 = (Q s q) — 1. In 
principle, such cut-off is proportional to the inverse of the coupling constant, but with 
an unknown proportionality constant which prevents us from making such connection 
with the value of a s in the definition of y. 



3 Results 

We have performed the evolution (j2J starting from the initial conditions indicated in 
Subsection 12.11 In Fig. ^ the evolution for both initial conditions GBW and MF, for 
9 given individual configurations with Q 2 = 0.01, 0.04, 0.1, 0.31, 1, 3.31, 10, 30 and 
100, for linear and logarithmic averages, and for the cut-off evolution 3 , are shown. 
The development of traveling waves [30] can be seen in all cases and for both initial 
conditions. For this reason, hereon we will present results only for MF, the results for 
GBW being in agreement. 



3 In the curves which result from this cut-off evolution, a discontinuity in the derivative induced by 
the cut-off is clearly visible. 
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Figure 1: BK evolution starting from GBW (upper-left plot) and MF (the other four 
plots) initial conditions, all for individual configurations with Q 2 sQ = 0.01, 0.04, 0.1, 
0.31, 1, 3.31, 10, 30 and 100 (thin lines left to right), for linear (A = 100, upper plots) 
and logarithmic averages (A/ = 100, middle-left plot, and A; = y, middle-right plot) 
using thick lines, and for cut-off evolution (thin lines, lower plot) for k = 0.05, 0.01 
and 0.002 (left to right). In this last plot, horizontal dotted lines indicate the values of 
k. Results are shown for y — 0, 2, 5 and 8 in black, red, green and blue respectively 
(sets of lines from left to right). 



7 



In these plots no dramatic effect of the averaging can be observed, the large-fc be- 
havior appearing very similar to that of the individual configurations. This fact could 
be expected from the power-like behavior of the solutions of BK [23,27,28,30] at large 
k. In order to see the details of the effect of the averaging we show in Fig. El zoomed 
versions of Fig. ^ in linear vertical scale. The pictures clearly show the mixing of 
individual configurations attained by the averaging procedure, substantiating for BK 
evolution the picture in [38] . It can also be seen that the effect of the averaging with a 
fixed dispersion is more evident in the initial condition and becomes less evident with 
increasing rapidities, in disagreement with present expectations from the sFKPP equa- 
tion as commented previously (see e.g. [51]), while for a dispersion linearly dependent 
on y, the effect becomes, in agreement with sFKPP predictions, more and more evident 
with increasing rapidities. 

In the next two Subsections we will examine in more detail two features of evolution, 
namely the approach to a universal scaling form and the evolution of the saturation 
scale with rapidity. 

3.1 Approach to scaling 

BK evolution is known [19,21,23,27,30] to approach asymptotically a scaling function, 

4>(y, k) — > <p(k/Q s (y)) for y — > oo, independent of the initial condition for evolution. On 

the other hand, attempts [37,38] to go beyond the mean field approximation underlying 

the BK equation predict a violation of the scaling behavior, leading to what has been 

called diffusive scaling [51]. Here we try to explore this aspect within our framework of 

modified BK evolution. For the function or the average of individual configurations, 

obtained by evolution up to a given y, we compute the saturation scale as the position of 

the maximum of the function h(k) = k 2 Vl<f)(k), in line with previous works [17, 19, 23] 4 . 

In standard BK evolution, this function shows a Gaussian-like shape which is preserved 

in the evolution, and the velocity of the evolution with y of this soliton-like form is 

given by the saturation scale. 

4 Note that this choice corresponds to a region where the function is large, 4> ~ b and not to a dilute 
region of the traveling wave. It is known [24,48,54] that with this definition, sub-leading terms [32] 
in the y-evolution of Q s are less noticeable. In any case, we do not attempt to study in detail the 
behavior of the saturation scale with rapidity but only to see whether clear differences among the 
different modifications of BK evolution appear or not. 
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Figure 2: Details of Fig. [T] (middle- left, upper plot here, and middle-right, lower plot 
here) in vertical linear scale. Results for y = 2 overlap with those for y = and are 
not shown for reasons of clarity. 



In Fig.|3]we show the results for evolution of function <fi when plotted versus kj Q s {y). 
Curves for an individual configuration with Q s q = 1, and for linear and logarithmic 
averages and for cut-off evolution, are shown. The violation of scaling due to averaging 
does not look dramatic, the scaling curves being relatively close and getting closer 
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for increasing rapidities. Only the very wide logarithmic averages seem to induce 
some visible violation of scaling, see below. On the other hand, scaling is achieved very 
quickly in the cut-off evolution (only for the largest cut-off and rapidity some departure 
is visible). This fact could be expected as it is known [24,27,30] that violations of scaling 
come from the tails of large k, which are cut in the cut-off evolution. 
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Figure 3: Scaling forms <j)(k/Q s (y)) for evolution starting from MF initial conditions, 
for y = 2, 5 and 8 (black, red and green lines respectively). Upper-left plot: results of 
an individual configuration with Q s0 = 1 (solid lines), of linear (A = 100, dashed lines) 
averaging and with cut-off (k = 0.01, dashed-dotted lines). Upper-right plot: idem 
but for logarithmic (A/ = 100) averaging. Lower-left plot: idem but for logarithmic 
(Ai = y) averaging. Lower-right plot: comparison of the scaling forms for different 
cut-offs k = 0.05 (solid lines), 0.01 (dashed lines) and 0.002 (dashed-dotted lines). In 
this last plot, horizontal dotted lines indicate the values of k. 



To analyze in more detail the effect of averaging on scaling, in Fig. ^zoomed versions 
in linear scale of the Fig.Elfor the logarithmic averages are shown. It is well known from 
standard BK evolution [19,21,23,27,30] that individual configurations tend towards 
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a universal scaling shape for the largest y, but violations of scaling of order ~ 10 % 
around Q s persist until large rapidities y ~ 10 [23]. The cut-off evolution shows some 
scaling violation, smaller than the standard BK results, and which tends to disappear 
with increasing k/Q s (y), see comments above. Besides, the averaged solutions for 
A^ = 100 show violations of scaling similar those of the individual solutions at the 
same y. But a dispersion increasing with rapidity, A; = y, produces larger violations 
of scaling than the case of standard BK evolution 5 . 

At first sight it looks odd that a smaller, varying dispersion creates larger violations 
of scaling than a larger, fixed dispersion. This is due to the fact that we sample a limited 
space of initial conditions, which for a very wide averaging gets immediately saturated, 
thus making the effect of the averaging less noticeable 6 . So we conclude that the 
effect on scaling of our cut-off evolution seems to contradict the results from sFKPP 
evolution (it should be noted that our cut-off is applied to an evolution performed on 
averages, not to the evolution of individual configurations), while the effect on scaling of 
our averaging procedure with a dispersion increasing with rapidity agrees qualitatively 
with the results from sFKPP. 

3.2 Rapidity evolution of the saturation scale 

To study the rapidity evolution of the saturation scale, we determine it as indicated in 
the previous Subsection. An error, corresponding to the distance between neighbors 
in the grid at the saturation scale, is assigned to every value of the saturation scale. 
In Fig. El we show the results for the evolution of the saturation scale with reduced 
rapidity y. In order to quantify eventual differences, we have performed a fit to an 
exponential form 7 Q 2 s {y) = ce dy in the region 3 < y < 9. 

The averaging procedure shows no sizable effect, a fact which contradicts the expec- 
tations from sFKPP but can be easily understood as BK evolution is known [19,21,23, 
5 From the comparison of the curves for rapidities y = 2, 5, 8 the violation of scaling induced by 
the logarithmic averaging with A; = y looks ~ 20 % for k/Q s {y) = 1, increasing to ~ 50 % for 

k/Q s (y) — 10, while the corresponding scaling violations for the individual solutions are smaller. 
6 This fact justifies our choice A; = y. 

7 The value of d expected [27,28,30] is d ~ 4.88; for evolution of individual configurations and for 
the average, we get a smaller value, d ~ 4.45, due to the neglected sub-leading terms [24,32,48,54] in 
the y-behavior of the saturation scale at these large but finite values of y. 
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Figure 4: Details of Fig. El (upper-right, upper plot here, and lower-left, lower plot 
here) in vertical linear scale, for values of k around Q s {y)- 



27,30] to lead asymptotically to a universal wave front moving with universal velocity 
(for 'supercritical' initial conditions, see the footnote in Subsection I2.1|) . Thus, asymp- 
totically every individual configuration (even those arising from initial conditions not 
only with different initial saturation scale but also with different functional shape) will 
move with the same velocity, and the average will be characterized by this common 
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Figure 5: Q 2 S versus y for evolution starting from MF initial conditions. In the upper 
and lower-left plots, the results for the evolution of individual configurations are shown 
for Q 2 s0 = 0.01, 1 and 100 in black, red and green respectively, while in blue the results 
for linear (A = 100, upper-left plot) and logarithmic (A; = 100, upper-right plot, and 
A; = y, lower-left plot) are presented. In the lower-right plot, results of the cut-off 
evolution for k = 0.002 (black), 0.01 (red) and 0.05 (green) are presented. In the plots, 
straight lines correspond to fits to Q 2 s {y) = ce dy in the region 3 < y < 9, with the 
values of d indicated in the plots. 



velocity. 

On the other hand, the cut-off evolution shows a strong influence on the ^-evolution 
of the saturation scale in agreement with expectations from sFKPP: the larger the value 
of the cut-off parameter k, the slower the evolution 8 . This corresponds to the fact that 
the evolution is driven [23, 27, 30] by the large- k regions, as discussed previously in 
Subsection 13. 11 



3 In the framework of sFKPP, behaviors proportional to (In k) 1 and k 1 are expected in the limit 



of weak [38,48,49] and strong [47,50] noise respectively. Our results show a mild behavior 



.-0.095 
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Other interesting quantity which we can examine 9 is not the saturation scale of the 
averaged solution, but the average saturation scale. In the language of [30-32], the 
evolution of the saturation scale with rapidity shows the speed of the wave front in 
its movement towards higher transverse momenta. Thus the dispersion of the satu- 
ration scale at a given rapidity will reflect the spread of the ensemble of wave fronts 
coming through evolution from a given ensemble of initial conditions. It is then a very 
good quantity to characterize the effect of the evolution on a possible transition from 
geometric to diffusive scaling [51]. 

We define the average of any observable (0(y)) by substituting (f>(k) by 0(y) in 
Eq. ©. In Fig. El we plot the average values (Ql(y)) and (In Q^(y)) and the corre- 
sponding dispersions, versus rapidity, for linear and logarithmic weight functions, for 
two different fixed values of A, A/, and also for A/ = y. In the case of fixed A, A; both 
the average Q 2 s {y) and its dispersion show the same dependence with y, so the width of 
the ensemble of wave fronts gets wider and wider with increasing rapidity. On the other 
hand, the dispersion of lnQg(y) stays constant (within our numerical accuracy) during 
evolution - a behavior opposite to the linear y-dependence found in sFKPP [38,47-50]. 
Besides, smaller initial widths result in smaller widths after evolution. We expect these 
characteristics to remain true for large enough rapidity even if the weight functions for 
averaging are more involved, or even for initial ensembles containing different shapes. 
This is due to the fact that BK evolution eventually leads to a universal shape which, 
once roughly built, will evolve with the characteristics we have discussed. On the other 
hand, for Ai = y the dispersion of In Q^(y) increases with increasing rapidity roughly 
linearly in agreement with expectations from sFKPP, though with a coefficient smaller 
than expected from the dispersion introduced in the weight function. 



4 Discussion 

Much theoretical effort has been directed recently to go beyond the standard JIMWLK 
framework to study the small-x structure of hadrons: the role of correlations [36,38,50, 
53,55,56], of discreteness in gluon emissions [38,57-59] and the relevance of Pomeron 
loops [37,39-46,60-64]. Numerical studies are just starting [25,47,48,59,65]. 
9 We thank Alex Kovner for drawing our attention on this point. 



14 



dashed! V(<Q S 4 > - <Q S 2 > 2 ) 
solid : <Q S 2 > 




dashed: V(<ln 2 Q s 2 >-<ln Q s 2 > 2 ) 
solid : < In Q 2 > 




y 

Figure 6: {Q 2 s {y)) (black solid lines in the plots on the left), (hxQ 2 s {y)) (black solid lines 
in the plots on the right), and their corresponding dispersions (red dashed lines in the 
plots on the left and right respectively) versus y, for linear (upper plots) and Gaussian 
(middle and lower plots) averaging with A, A; = 100 (thick lines), and A = 1, A/ = 1 
and Ai = y (thin lines in the upper, middle and lower plots respectively), for evolution 
starting from MF, 6 = 1, initial conditions. The difference between thick and thin red 
dashed lines in the upper-right plot, and between thick and thin black solid lines in 
the middle- and lower-right plots, is numerically very small and hardly visible. 
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In this paper we contribute to this subject through an attempt to quantify numer- 
ically the relevance of two of these aspects, namely correlations through averaging [53] 
and discreteness of gluon emissions [38], on standard BK evolution. The averaging of 
individual configurations results in scaling violations, but their size around the satura- 
tion scale is not much larger than the violations already present in standard evolution 
at a large but finite rapidity. On the other hand, the evolution of the saturation scale 
with rapidity, and the large-A; behavior of the solutions, seem to be hardly affected. 
Finally, the dispersion in the saturation scale of the individual solutions increases with 
rapidity at least as fast as the corresponding average, indicating that the width of the 
ensemble of evolved solutions gets wider and wider with increasing rapidity. 

The implementation of discreteness in the evolution through a cut-off does not 
induce any further scaling violation. On the other hand, the saturation scale evolves 
slower with rapidity with increasing value of the cut-off, and the solutions at large k 
become steeper than in the standard BK evolution. 

When compared with the expectations from the sFKPP equation [39,42,47-50], we 
conclude than none of the proposed modifications is able to reproduce all the features 
found there. The implementation of a cut-off slows down the evolution in accordance 
with expectations, but results in a better fulfillment of scaling, at variance with sFKPP. 
The averaging procedure, when implemented through a dispersion increasing with ra- 
pidity, leads to 'diffusive' wave fronts [51], but the rapidity evolution of the saturation 
scale is not affected, contrary to expectations. 

In view of the differences among the effects of different modifications on stan- 
dard BK evolution and with the results of sFKPP, of the fact that the experimentally 
accessible regions of rapidity are pre-asymptotic and dominated by the initial condi- 
tions [24,34,35], and of the additional uncertainties on the form of the averaging weight 
function and/or size of the cut-off, any claim on phenomeno logical implications of our 
study looks premature. Besides, our implementation of these modifications of BK evo- 
lution is simplistic, so their results could only be taken as hints of the possible effects 
of a proper application of such modifications. Nevertheless, we find their effects large 
and diverse enough to justify the further investigations, both on the analytical and on 
the numerical side, currently under development to rigorously develop and implement 
all these new ideas. Work in this direction is in progress. 
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